! (C) Copyright 2005- ECMWF.
!
! This software is licensed under the terms of the Apache Licence Version 2.0
! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
! In applying this licence, ECMWF does not waive the privileges and immunities
! granted to it by virtue of its status as an intergovernmental organisation
! nor does it submit to any jurisdiction.

PROGRAM CREATE_INIT_CLIM
USE NETCDF
USE GET_POINTS_MOD

INTEGER :: NCID_IN,NCID_OUT,NCID,VARID,DIMID

INTEGER :: NP_TOT,NPOINTS,NMONTH,NLEVS
INTEGER     ,DIMENSION(:),ALLOCATABLE :: POINTS 
REAL(KIND=8),DIMENSION(:),ALLOCATABLE :: LAT,LON
CHARACTER(LEN=500) :: INIT_OUT,INIT_IN,CLIM_IN,CLIM_OUT,INFO_FILE,NAM_FILE,FILE_POINTS
INTEGER :: NP_LAND_OUT,NP_LAKE_OUT,NP_LAKE_FRAC_OUT,NP_OCEAN_OUT
REAL(KIND=8)       :: LDEPTH_IN
LOGICAL :: LINIT,LCLIM

NAMELIST /NAMFP/INFO_FILE
NAMELIST /NAMCC/LDEPTH_IN,NMONTH,NLEVS,LINIT,LCLIM,INIT_OUT,INIT_IN,CLIM_OUT,CLIM_IN,FILE_POINTS

NAM_FILE='input.nam'
! Default Values 
LDEPTH_IN=30. ! DEFAULT VALUE 
NMONTH=12     ! number of months 
NLEVS=4       ! DEFAULT VALUE (10-layer experimental)
LINIT=.True.
LCLIM=.True. 
INIT_OUT='surfinit.nc'
INIT_IN='surfinit_all.nc'
CLIM_OUT='surfclim.nc'
CLIM_IN='surfclim_all.nc'
FILE_POINTS='NONE'

OPEN(99,FILE=TRIM(NAM_FILE),STATUS='old')
READ(99,NAMFP)
READ(99,NAMCC)
CLOSE(99)

IF (FILE_POINTS == 'NONE') THEN
    CALL GET_POINTS(TRIM(INFO_FILE),TRIM(NAM_FILE),NP_TOT,NPOINTS,POINTS,&
                & NP_LAND_OUT,NP_LAKE_OUT,NP_LAKE_FRAC_OUT,NP_OCEAN_OUT)
	    
    PRINT*,'NP_TOT',NP_TOT
    PRINT*,'NPOINTS',NPOINTS
    PRINT*,'NP_LAND',NP_LAND_OUT
    PRINT*,'NP_LAKE',NP_LAKE_OUT
    PRINT*,'NP_LAKE_FRAC',NP_LAKE_FRAC_OUT
    PRINT*,'NP_OCEAN',NP_OCEAN_OUT
ELSE
    CALL NCERROR( NF90_OPEN(FILE_POINTS,NF90_NOWRITE,NCID) )
    CALL NCERROR( NF90_INQ_DIMID(NCID,'x',DIMID) )
    CALL NCERROR( NF90_INQUIRE_DIMENSION(NCID,DIMID,LEN=NPOINTS) )
    
    ALLOCATE(POINTS(NPOINTS))
    
    CALL NCERROR( NF90_INQ_VARID(NCID,'x',VARID) )
    CALL NCERROR( NF90_GET_VAR(NCID,VARID,POINTS) )
    
    CALL NCERROR( NF90_CLOSE(NCID) )

ENDIF 

IF (LINIT) THEN 
    !=============================================
    ! INIT 
    ! OPEN FILES 
    CALL NCERROR( NF90_OPEN(INIT_IN,NF90_NOWRITE,NCID_IN) )
    CALL NCERROR( NF90_OPEN(INIT_OUT,NF90_WRITE,NCID_OUT) )
    
    IF (FILE_POINTS /= 'NONE') THEN
        CALL NCERROR( NF90_INQ_DIMID(NCID_IN,'x',DIMID) )
        CALL NCERROR( NF90_INQUIRE_DIMENSION(NCID,DIMID,LEN=NP_TOT) )
    ENDIF 
    
    PRINT*,''
    PRINT*,'**********************************************'
    PRINT*,'TREATING INIT'
    PRINT*,'READING FROM->',TRIM(INIT_IN)
    PRINT*,'WRITING TO  ->',TRIM(INIT_OUT)

    CALL WRITE2CDF(NCID_IN,NCID_OUT,NP_TOT,NPOINTS,POINTS,NMONTH,NLEVS)

    ! CLOSE FILES 
    CALL NCERROR( NF90_CLOSE(NCID_IN))
    CALL NCERROR( NF90_CLOSE(NCID_OUT))
ENDIF

IF (LCLIM) THEN
    !=============================================
    ! CLIM
    ! OPEN FILES 
    CALL NCERROR( NF90_OPEN(CLIM_IN,NF90_NOWRITE,NCID_IN) )
    CALL NCERROR( NF90_OPEN(CLIM_OUT,NF90_WRITE,NCID_OUT) )
    
    IF (FILE_POINTS /= 'NONE') THEN
        CALL NCERROR( NF90_INQ_DIMID(NCID_IN,'x',DIMID) )
        CALL NCERROR( NF90_INQUIRE_DIMENSION(NCID,DIMID,LEN=NP_TOT) )
    ENDIF 
    
    PRINT*,''
    PRINT*,'**********************************************'
    PRINT*,'TREATING CLIM'
    PRINT*,'READING FROM->',TRIM(CLIM_IN)
    PRINT*,'WRITING TO  ->',TRIM(CLIM_OUT)

    CALL WRITE2CDF(NCID_IN,NCID_OUT,NP_TOT,NPOINTS,POINTS,NMONTH,NLEVS)

    ! CLOSE FILES 
    CALL NCERROR( NF90_CLOSE(NCID_IN))
    CALL NCERROR( NF90_CLOSE(NCID_OUT))
ENDIF 

CONTAINS
SUBROUTINE WRITE2CDF(NCID_IN,NCID_OUT,NP_TOT,NPOINTS,POINTS,NMONTH,NLEVS)
USE NETCDF
IMPLICIT NONE

INTEGER,INTENT(IN) :: NCID_IN,NCID_OUT,NP_TOT,NPOINTS,NMONTH,NLEVS
INTEGER,DIMENSION(NPOINTS),INTENT(IN) :: POINTS 

INTEGER  :: NVARS,IVAR,JP,VARID,STATUS,VARID1
CHARACTER(LEN=50) :: VNAME,VNAME1
REAL(KIND=8) :: DATA_IN(NP_TOT),DATA_IN1(NP_TOT),DATA_IN2(NP_TOT,NLEVS),DATA_INM(NP_TOT,NMONTH)
REAL(KIND=8) :: DATA_OUT(NPOINTS),DATA_OUT2(NPOINTS,NLEVS),DATA_OUTM(NPOINTS,NMONTH),HLICEMIN
REAL :: DATA_VALUE

HLICEMIN=0.001D0  ! ALL VALUES BELLOW THIS SET TO zero!
CALL NCERROR( NF90_INQUIRE(NCID_OUT,nVariables=NVARS) )
PRINT*,'WRITE2CDF'
PRINT*,'TREATING NVARIABLES=',NVARS

DO IVAR=1,NVARS
   
   CALL NCERROR( NF90_INQUIRE_VARIABLE(NCID_OUT,IVAR,name=VNAME) )
    VARID=-1
   SELECT CASE (VNAME)
    CASE ('x')
      DO JP=1,NP_TOT
        DATA_IN(JP)=REAL(JP)
      ENDDO
    CASE ('lat')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'lat',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('lon')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'lon',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('SWE')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SD',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
      DATA_IN(:)=DATA_IN(:)*1000. 
      PRINT*,'SD (m) RESCALLED TO kg.m-2'
    CASE ('snowdens')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'RSN',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('SAlbedo')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ASN',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('SnowT')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'TSN',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('CanopInt')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SRC',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('AvgSurfT')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SKT',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('TLICE')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LICT',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LICD',VARID1) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID1,DATA_IN1,(/1/),(/NP_TOT/)) )
      DO JP=1,NP_TOT
        DATA_IN(JP)=DMIN1(DATA_IN(JP),273.15D0)
        IF ( DATA_IN1(JP) < HLICEMIN ) THEN
            DATA_IN(JP) = 273.15D0
        ENDIF
      ENDDO
    CASE ('TLMNW')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LTLT',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
      DO JP=1,NP_TOT
        DATA_IN(JP)=DMAX1(DATA_IN(JP),273.15D0)
      ENDDO
    CASE ('TLWML')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LMLT',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
      DO JP=1,NP_TOT
        DATA_IN(JP)=DMAX1(DATA_IN(JP),273.15D0)
      ENDDO
    CASE ('TLBOT')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LBLT',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
      DO JP=1,NP_TOT
        DATA_IN(JP)=DMAX1(DATA_IN(JP),273.15D0)
      ENDDO
    CASE ('TLSF')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LSHF',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
      DO JP=1,NP_TOT
        DATA_IN(JP)=DMAX1(DATA_IN(JP),0.5D0)
      ENDDO
    CASE ('HLICE')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LICD',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
      DO JP=1,NP_TOT
        DATA_IN(JP)=DMAX1(DATA_IN(JP),0.0D0)
        IF ( DATA_IN(JP) < HLICEMIN ) THEN 
            DATA_IN(JP)=0.D0
        ENDIF
      ENDDO
    CASE ('HLML')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LMLD',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
      DO JP=1,NP_TOT
        DATA_IN(JP)=DMIN1(DATA_IN(JP),50.0D0)
      ENDDO
    CASE ('SoilTemp')
      IF (NLEVS == 4) THEN
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,1),(/1,1/),(/NP_TOT,1/)),'rr1' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL2',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,2),(/1,1/),(/NP_TOT,1/)),'rr2' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL3',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,3),(/1,1/),(/NP_TOT,1/)),'rr3' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,4),(/1,1/),(/NP_TOT,1/)),'rr4' )
      ELSEIF (NLEVS == 10) THEN
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,1),(/1,1/),(/NP_TOT,1/)),'rr1' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,2),(/1,1/),(/NP_TOT,1/)),'rr1' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,3),(/1,1/),(/NP_TOT,1/)),'rr1' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL2',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,4),(/1,1/),(/NP_TOT,1/)),'rr2' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL2',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,5),(/1,1/),(/NP_TOT,1/)),'rr2' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL3',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,6),(/1,1/),(/NP_TOT,1/)),'rr3' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL3',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,7),(/1,1/),(/NP_TOT,1/)),'rr3' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,8),(/1,1/),(/NP_TOT,1/)),'rr4' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,9),(/1,1/),(/NP_TOT,1/)),'rr4' )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'STL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,10),(/1,1/),(/NP_TOT,1/)),'rr4' )      
      ELSE
        PRINT*,'Vertical discretisation not supported for NLEVS=',NLEVS
      ENDIF
    CASE ('icetemp') 
      IF (NLEVS == 4) THEN
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,1),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL2',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,2),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL3',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,3),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,4),(/1,1/),(/NP_TOT,1/)) )
      ELSEIF (NLEVS == 10) THEN
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,1),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,2),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,3),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL2',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,4),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL2',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,5),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL3',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,6),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL3',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,7),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,8),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,9),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'ISTL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,10),(/1,1/),(/NP_TOT,1/)) )
      ELSE
        PRINT*,'Vertical discretisation not supported for NLEVS=',NLEVS
      ENDIF
    CASE ('SoilMoist')
      IF (NLEVS == 4) THEN
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,1),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL2',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,2),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL3',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,3),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,4),(/1,1/),(/NP_TOT,1/)) )
      ELSEIF (NLEVS == 10) THEN
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,1),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,2),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL1',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,3),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL2',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,4),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL2',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,5),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL3',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,6),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL3',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,7),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,8),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,9),(/1,1/),(/NP_TOT,1/)) )
        CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SWVL4',VARID) )
        CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN2(:,10),(/1,1/),(/NP_TOT,1/)) )
      ELSE
        PRINT*,'Vertical discretisation not supported for NLEVS=',NLEVS
      ENDIF
    CASE ('Mask')
      DATA_IN(:)=1.
    CASE ('Malbedo')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'AL',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_INM,(/1,1/),(/NP_TOT,NMONTH/)) )
    CASE ('Mlail')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LAIL',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_INM,(/1,1/),(/NP_TOT,NMONTH/)) )
    CASE ('Mlaih')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LAIH',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_INM,(/1,1/),(/NP_TOT,NMONTH/)) )
    CASE ('z0m')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SR',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('lz0h')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LSRH',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('landsea')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'LSM',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('geopot')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'Z',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('cvl')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'CVL',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('cvh')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'CVH',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('tvl')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'TVL',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('tvh')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'TVH',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) ) 
    CASE ('sdor')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SDOR',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('sst')
!       CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SSTK',VARID) )
!       CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
        DATA_IN(:)=280.
    CASE ('seaice')
!       CALL NCERROR( NF90_INQ_VARID(NCID_IN,'CI',VARID) )
!       CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
        DATA_IN(:)=0.
    CASE ('sotype')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'SLT',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
      DO JP=1,NP_TOT
        IF (DATA_IN(JP) == 0 ) DATA_IN(JP)=1.
      ENDDO
    CASE ('LDEPTH')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'DL',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE ('CLAKE')
      CALL NCERROR( NF90_INQ_VARID(NCID_IN,'CL',VARID) )
      CALL NCERROR( NF90_GET_VAR(NCID_IN,VARID,DATA_IN,(/1/),(/NP_TOT/)) )
    CASE DEFAULT 
      PRINT*, 'VARIABLE NOT FOUND !!!'
      STOP 9
   END SELECT 
   SELECT CASE (VNAME)
     CASE ('icetemp','SoilMoist','SoilTemp') 
       DO JP=1,NPOINTS
         DATA_OUT2(JP,1:NLEVS)=DATA_IN2(POINTS(JP),1:NLEVS)
       ENDDO 
       CALL NCERROR( NF90_PUT_VAR(NCID_OUT,IVAR,DATA_OUT2(1:NPOINTS,1:NLEVS),(/1,1/),(/NPOINTS,NLEVS/) ) )
     CASE ('Malbedo','Mlail','Mlaih')
       DO JP=1,NPOINTS
         DATA_OUTM(JP,1:NMONTH)=DATA_INM(POINTS(JP),1:NMONTH)
       ENDDO 
       CALL NCERROR( NF90_PUT_VAR(NCID_OUT,IVAR,DATA_OUTM(1:NPOINTS,1:NMONTH),(/1,1/),(/NPOINTS,NMONTH/) ) )
     CASE DEFAULT
       DO JP=1,NPOINTS
         DATA_OUT(JP)=DATA_IN(POINTS(JP))
       ENDDO 
       CALL NCERROR( NF90_PUT_VAR(NCID_OUT,IVAR,DATA_OUT(1:NPOINTS),(/1/),(/NPOINTS/) ) )
       
   END SELECT 
   
   STATUS=NF90_INQUIRE_VARIABLE(NCID_IN,VARID,name=VNAME1) 
   IF( STATUS /= 0 ) THEN
      PRINT*,'DEFAULT->',TRIM(VNAME)
      PRINT*,TRIM(VNAME),'==',DATA_OUT(1)
   ELSE
    PRINT*,TRIM(VNAME1),'->',TRIM(VNAME)
   ENDIF

ENDDO


END SUBROUTINE WRITE2CDF

END PROGRAM CREATE_INIT_CLIM

